#!/usr/bin/python
#-*-coding: utf-8 -*-

'''
Created on 11 mai 2012

@author: puiseux
'''
import sys
import numpy as np
from path import Path
import pickle
from exceptions import *
import datetime
import scipy as sp
# from scipy import interpolate, sum, prod, ones
from random import choice
from PyQt4 import QtGui, QtCore
import math

DEBOG = True
VERBOSITY_LEVEL = 2
DEBOG_OUTPUT = sys.stderr 
DEBOG_OUTPUT = sys.stdout
output = DEBOG_OUTPUT

def rotate(points, alfa, centre):
    '''
     Retourne une COPIE de points rotationnée(!).  
    points est supposé stocké par ligne (shape=(n,2)), chaque point est de shape (1,2), il les faudrait en colonne (shape=(2,1)) pour faire le produit matriciel. 
    Donc on transpose tout et on ecrit Xi' = C' + (Xi'-C')*A' au lieu de Xi = C + A*(Xi-C), pour i= 0, 1,...
    ''' 
    Ct = np.asarray(centre).reshape((1, 2))
    cosa, sina = np.cos(alfa), np.sin(alfa)
    At = np.matrix([[cosa, -sina],[sina, cosa]]).transpose()
    Xt = points-Ct
    Xt =  Xt*At + Ct
#     trace('', Xt.shape)
    return np.asarray(Xt)


def symetrieAxe(points, axe, centre=None):
    '''
    homothétie de centre 'centre': X = C + ech*(X-C) = (1-ech)*C + ech*X
    points est modifié in situ.
    '''
#     trace(None, points.shape, axe, centre)
    points[:,axe] *= -1
    points[:,axe] += 2*centre
    return points

def hardScale(points, echelle, centre=None, translate=False):
    '''
    homothétie de centre 'centre': X = C + ech*(X-C) = (1-ech)*C + ech*X
    points est modifié in situ.
    '''
    if centre is None : 
        centre = [0,0]
    centre = np.asarray(centre).reshape((1,2))
    points *= [echelle[0], echelle[1]] 
    points += [1.0-echelle[0], 1.0-echelle[1]]*centre
    return points
 
def baryCentre(points, masses=None):
    '''Retourne le barycentre du nuage de points 'points' affectés des masses 'masses' 2d'''
    points = np.asarray(points)
#     trace('', points.shape, len(points))
    N = len(points)
    if masses is None : 
        masses = np.ones((N,))/N
    else : 
        masses = masses/np.sum(masses)
    X, Y = np.prod([masses, points[:,0]], axis=0), np.prod([masses, points[:,1]], axis=0) 
    cg = [[np.sum(X), np.sum(Y)],]
    return np.asarray(cg)
 

def p2t(qpoint):
    """ retourne les coord d'un QPointF sous forme de tuple simple"""
    return (qpoint.x(), qpoint.y())

def p2a(qpoint):
    """ retourne les coord d'un QPointF sous forme de np.ndarray((1,2))"""
    return np.asarray((qpoint.x(), qpoint.y())).reshape((1,2))

def p2s(qpoint):
    """ retourne les coord d'un QPointF sous forme de chaine"""
    return '(%.2f, %.2f)'%(qpoint.x(), qpoint.y())

def dist2(p1, p2, n=2):
    '''retourne le carré de la distance de p1 à p2 en norme n=2'''
    if isinstance(p1, QtCore.QPointF) : x1, y1 = p1.x(), p1.y()
    else : x1, y1 = p1[0], p1[1]
    if isinstance(p2, QtCore.QPointF) : x2, y2 = p2.x(), p2.y()
    else : x2, y2 = p2[0], p2[1]
    return (x2-x1)**2+(y2-y1)**2
    
def dist(p1, p2, n=2):
    '''retourne la distance de p1 à p2 en norme n=2'''
    return math.sqrt(dist2(p1, p2))


def absCurv(points):
    """Abscisse curviligne des points de points"""
    T = np.zeros(len(points))
    for k in range(len(points)-1) :
        T[k+1] = T[k] + dist(points[k],points[k+1])
    return T


def splineCubique(points, weights=None, bbox=((0,-1),(0,-1)), eps=1.0e-2, degre=(2,3)):
    """
    Ajustement d'une spline cubiques paramétriques au polygone points
    sx(t), sy(t) sont deux splines.
    Retourne:
    -------- 
    T = les valeurs du paramètre t, abscisse curviligne NORMALISEE, entre 0 et 1.
    sy, sy = les deux splines
    TODO: parametres de précision et d'ajustement  
    - ajust : proportionnel à la surface (au lieu de T[-1]*T[-1])
    - s : proportionnel à deltax et deltay suivant le cas (au lieu de T[-1])
    """
    N = len(points)
    if len(points) <3 :
        raise PyGlideSplineError('Pas assez de points (%d<3), lissage impossible'%N)
    T = absCurv(points)
    T /= T[-1]
    #print whoami(), 'len(points), T[-1]', len(points), T[-1]
    ajust = eps*T[-1]*T[-1]
    for k, weight in enumerate(weights) : 
        weights[k] = weight if weight>=0 else len(weights)
    weights *= N/sum(weights)#Charge totale = N
    bbx = [T[bbox[0][0]], T[bbox[0][1]]]
    bby = [T[bbox[1][0]], T[bbox[1][1]]]
    try : 
        sx = sp.interpolate.UnivariateSpline(T, points[:,0], w=weights, bbox=bbx, k=degre[0], s=ajust*T[-1])#s=la précision de l'ajustement s=0 <=> interpolation
        sy = sp.interpolate.UnivariateSpline(T, points[:,1], w=weights, bbox=bby, k=degre[1], s=ajust*T[-1])
    except : 
        raise PyGlideSplineError("%s : impossible de calculer la spline"%whoami())
    return T, sx, sy

def aire(a,b,c):
    """
    Aire du triangle abc dans l'espace.
    C'est la moitié de la norme du produit vectoriel ab vect ac
    """
    u, v = b-a, c-a
    r = u[2]*v[0]-u[0]*v[2]
    s = u[0]*v[1]-u[1]*v[0]
    t = u[1]*v[2]-u[2]*v[1]
    return 0.5*np.sqrt(r*r+s*s+t*t)

def segmentPlusProche(points, P):
        u'''retourne le numeros i du segment points[i, i+1] le plus proche de P. 
        1- On cherche le point le plus proche
        2- on regarde de quel coté on est
        voir sur http://www.exaflop.org/docs/cgafaq/cga1.html#Subject%201.02:%20How%20do%20I%20find%20the%20distance%20from%20a%20point%20to%20a%20line?
        '''
        '''Subject 1.02: How do I find the distance from a point to a line?

Let the point be C (Cx,Cy) and the line be AB (Ax,Ay) to (Bx,By).    The length of the line segment AB is L:

        L= sqrt( (Bx-Ax)^2 + (By-Ay)^2 ) .

Let P be the point of perpendicular projection of C onto AB. Let r be a parameter to indicate P's location along the line containing AB, with the following meaning:

          r=0      P = A
          r=1      P = B
          r<0      P is on the backward extension of AB
          r>1      P is on the forward extension of AB
          0<r<1    P is interior to AB

Compute r with this:

            (Ay-Cy)(Ay-By)-(Ax-Cx)(Bx-Ax)
        r = -----------------------------
                        L^2

The point P can then be found:

        Px = Ax + r(Bx-Ax)
        Py = Ay + r(By-Ay)

And the distance from A to P = r*L.

Use another parameter s to indicate the location along PC, with the following meaning:

           s<0      C is left of AB
           s>0      C is right of AB
           s=0      C is on AB

Compute s as follows:

            (Ay-Cy)(Bx-Ax)-(Ax-Cx)(By-Ay)
        s = -----------------------------
                        L^2

Then the distance from C to P = s*L.
        '''

        if isinstance(P, (QtCore.QPointF, QtCore.QPoint,)):
            X , Y  = P.x(), P.y()# Point argument
        else :
            X, Y = P[0], P[1]
        distances = [dist(point, (X,Y)) for point in points]
#        trace(None, 'P', X,Y)
#        trace(None, 'distances', distances)
        i0 = np.argmin(distances)
#        trace(None, i0)
        im, ip = (i0-1)%len(points), (i0+1)%len(points)
        #distance de P au segment [i,im]
#        trace(None, 'point le plus proche =', points[i0])
        x0, y0 = p0 = points[i0]# Point le plus proche
        xm, ym = pm = points[im]
        llm = dist2(p0,pm)
        lm = dist(p0,pm)
#        trace(None, 'point le plus proche =(%.2f,%.2f)'%(p0[0], p0[1]))
        rm = ((x0-X)*(x0-xm)+(y0-Y)*(y0-ym))/llm#produit scalaire
        dm = abs((x0-X)*(y0-ym)-(y0-Y)*(x0-xm))/lm#produit vectoriel, #distance P au segment
        xp,yp = pp = points[ip]
        llp = dist2(p0,pp)
        rp = ((x0-X)*(x0-xp)+(y0-Y)*(y0-yp))/llp#produit scalaire
        lp = dist(p0,pp)
        dp = abs((x0-X)*(y0-yp)-(y0-Y)*(x0-xp))/lp#produit vectoriel
#        print '(rm, rp)=(%g, %g) et (dm, dp)=(%g, %g)'%(rm, rp, dm, dp)
        if 0 <= rm <= 1 and 0 <= rp <= 1 : 
#            print 'in, in',
            if dm < dp :
#                print '=> -' 
                return im
            else: 
#                print '=> +' 
                return i0 #if i0!=0 else im
        elif 0 <= rm <= 1 and not (0 <= rp <=1) :
#            print 'in, out => -'
            return im
        elif 0 <= rp <= 1 and not (0 <= rm <=1) :
#            print 'out, in => +'
            return i0
        elif dm < dp : 
#            print 'out, out => -'
            return i0 
        else : 
#            print 'out, out => +'
            return im

def encombrement(piece):
    """
    piece doit être un np.ndarray
    retourne le paralelepipede d'encombrement d'un nuage de points 3d de la forme (N,3)
    """
    if isinstance(piece, (np.ndarray,)):
        points = piece.view()
        points.shape = points.size/3,3
    #    x, y, z = X,Y,Z = points[0,:]
        Max, Min = np.max, np.min #ca reste local
        M = np.asarray([Max(points[:,0]), Max(points[:,1]), Max(points[:,2])])
        m = np.asarray([Min(points[:,0]), Min(points[:,1]), Min(points[:,2])])
            
        return(m, M)
    elif isinstance(piece, (list, tuple, )) : 
        #liste de points 2d [(0.0,0.0) , (1.0,0.0),...]
        xmax = xmin = piece[0][0]
        ymax = ymin = piece[0][1]
        for point in piece : 
            xmax, xmin = max(xmax, point[0]), min(xmin, point[0])
            ymax, ymin = max(ymax, point[1]), min(ymin, point[1])
        return (xmin, ymin), (xmax, ymax)
    
    
def maintenant():
    u'''La date et l'heure formatées "human readable"'''
    return datetime.datetime.now().strftime("%A, %d. %B %Y %I:%M%p")

def doublons(points, voisinages=[], eps=1.0e-8, sym=True):
    """
    Verifier qu'il n'y a pas de point double.
    Paramètres : 
    ----------
    - points : ndarray(nbp,3) nbp points 3d. nbp peut être une 
        shape (nl, nc, 3) et dans ce cas nbp = nl*nc
    - voisinages ndarray(shape=(np,1+nv), dtype=int) np numéros de points, 
        faisant référence au tableau 'points', chacun des np points à 
        vérifier possède nv voisins. voisinage est donc de la forme
        [
         [k_0,  v_0_1,  v_0_2, ...,  v_0_nv],
         [k_1,  v_1_1,  v_1_2, ...,  v_1_nv],
         ...
         [k_np, v_np_1, v_np_2, ..., v_np_nv]
        ]
        les valeurs k_i, v_i_j sont des entiers entre 0 et nbp-1
    - eps : réel. deux points p,q sont considérée comme confondus si dist(p,q)<eps
    - sym : booleen. True si i voisin de j <==> j voisin de i. Divise par deux 
        le travail.
    
    Retourne :
    --------
    une liste de doublons (d_0,dd_0), ... (d_n, dd_n)
    """
    def cond(v,k):
        if sym : return v > k
        else : return v != k
        
    qoints = points.view()
    qoints.shape = (qoints.size/3,3)
    nbp = qoints.shape[0]
    if voisinages != [] : 
        voisinages = np.asarray(voisinages)
#        debog(whoami(), voisinages.shape)
        nbp = voisinages.shape[-1]
        voisinages.shape = (-1,nbp)
    else :
        if sym :  
            voisinages = np.asarray([[i]+range(i,nbp) for i in range(nbp)])
        else :    
            voisinages = np.asarray([[i]+range(nbp) for i in range(nbp)])

    doublons = []
    for voisinage in voisinages:
        k, voisins = voisinage[0], voisinage[1:]
        point = qoints[k]
        for voisin in voisins :
            if cond(voisin, k) : 
                if np.linalg.norm(point - qoints[voisin]) < eps :
                    doublons.append((k, voisin))
    return doublons


def sexplore(objet):
    '''retourne une chaine de caracteres avec les attributs de objet et leur valeur'''
    s = []
    dexp = explore(objet)
    for key in sorted(dexp.keys()) :
        s.append('%20s : %s'%(key, dexp[key]))
    return '\n'.join(s)

def explore(objet):
    '''retourne un dictionnaire avec les attributs de objet et leur valeur'''
    d = {}
    for key, value in objet.__dict__.iteritems() :
        if type(value) in (list, tuple) :
            try : d[key] = "[%s]"%(type(value[0]).__name__)
            except IndexError: d[key] = ( "[]")
        else : 
            d[key] ='%s'%(type(value).__name__)
    return d

def debog(*args, **kargs):
    if not DEBOG : return
    sys.stdout.flush()
    if VERBOSITY_LEVEL>1 :
        print>>output, "&BEGIN DEBOG"
    for arg in args : 
        print>>output, arg
    for key, value in kargs : 
        print>>output,key,' : ', value
    if VERBOSITY_LEVEL>1 :
        print>>output, "&END DEBOG"
            
def trace(objet=None, *args, **kargs):
    u'''L'instruction qui suit, dans une méthode de LaClasse 
        écrit sur sys.stdout, la provenance de l'appel (LaClasse::une_methode) puis une_variable, une_autre_variable
        (autant de variable que l'on veut)
    >>> trace(self, 'machin, truc', une_variable, une_autre_variable)
    '''
#    if not DEBOG : return
    if objet is None: 
        msg = phrase()
    else:
        try : msg = objet.__class__.__name__
        except AttributeError : msg=''
    lmsg = [msg+'::'+sys._getframe(1).f_code.co_name]
    for arg in args : 
        lmsg.append(str(arg))
    for key, value in kargs : 
        lmsg.append(key+' : '+ str(value))
        
    print>>output, ' ; '.join(lmsg)
    
def alert(objet=None, *args, **kargs):
    u'''Comme trace, mais sur sys.stderr (en rouge par défaut)'''
    if objet is None: 
        msg = ''
    else:
        try : msg = objet.__class__.__name__+'.'
        except AttributeError : msg=''
    lmsg = [msg+'::'+sys._getframe(1).f_code.co_name]
    for arg in args : 
        lmsg.append(str(arg))
    for key, value in kargs : 
        lmsg.append(key+' : '+ str(value))
        
    print>>sys.stderr, ' ; '.join(lmsg)
    
#alert = trace#TODO: faire un QMessageBox
      
def whoami(objet=None):
    if objet is None: 
        msg = ''
    else:
        try : msg = objet.__class__.__name__+'.'
        except AttributeError : msg=''
        
    return '%s::'%(msg + sys._getframe(1).f_code.co_name)#, id(object)

def toDict(cles, valeurs):
    """
    - cles = "cle1 cle2 cle3 ..." 
    - valeurs = "val1 val2 val3...", les valeurs sont des entiers ou des reels 
    retourne un dictionnaire cle,valeurs
    """
    d = {}
    for key, value in zip(cles.split(),valeurs.split()) :
        try: w = int(value)
        except ValueError : w = float(value)
        d[key.lower()] = w
    return d

def findAll(tag, lines, first=0, last=-1):
    """
    Retourne une liste triée des numeros de ligne de 
    TOUTES les occurences de tag dans lines[first:last]
    """
    #
    n = first -1
    n0 = _findRel (tag, lines, n+1, last)
    N = []
    while n0 >= 0 : 
        n = n0 + n + 1
        N.append(n)
        n0 = _findRel(tag, lines, n+1, last)
    return tuple(N)

def findAllLines(tag, lines, first=0, last=-1):
    '''
    Comme findAll(), mais il faut que la ligne complete (nettoyée) soit égale à tag, pas seulement une partie de la ligne.
    Par exemple : line = 'TOTO_EST_CONTENT' ne match pas avec tag='TOTO'
    '''
    liste = findAll(tag, lines, first, last)
    newlist = []
    for n in liste : 
        if lines[n].strip() == tag : 
            newlist.append(n)
    return tuple(newlist)
 
def _findRel(tag, lines, first=0, last=-1):
    """
    Cherche 'tag' dans lines[first:last] et retourne 
        * le numero de ligne i-first de la premiere occurence trouvee, 
          c'est à dire le numéro de ligne RELATIF : dans self.lines[0:]
        * None si le tag n'est pas trouvé
    """
    found = find(tag, lines, first, last)
    if found is None : return None
    else : return found - first
    
def find0(tag, lines, first=0, last=-1):
    """
    Cherche 'tag' dans lines[first:last] et retourne 
        * le numero de ligne i-first de la premiere occurence trouvee, 
          c'est à dire le numéro de ligne RELATIF : dans self.lines[0:]
        * None si le tag n'est pas trouvé
    """
    found = _findRel(tag, lines, first=first, last=last)
    if found is not None : return  found + first
    else : return None 


def find(tag, lines, first=0, last=-1):
    """
    Cherche 'tag' dans lines[first:last] et retourne 
        * le numero de ligne i de la premiere occurence trouvee, 
          c'est à dire le numéro de ligne ABSOLU : dans self.lines[0:]
        * None si le tag n'est pas trouvé
    """
    if last < 0 : last = len(lines)+last
    if first is None : first = 0 
    elif first < 0 : first = len(lines) + first
    if not tag : return None 
    i = first-1
    while i < last:
        i += 1 
        try : 
            if lines[i].find(tag) >= 0 :
                return i
        except IndexError : 
            return None
    return None

def rreplace(chaine, sub, replace):
    '''
    Remplace récursivement dans chaine, toutes les occurences de sub  par replace
    >>> rreplace("aaaaab","aa","a")
    ... "ab"
    '''
    while 1:
        oc = chaine
        chaine = oc.replace(sub, replace)
        if oc==chaine : break
    return chaine
    
def kompare(f1,f2):

    """
    Comparaison des variables de deux fichiers cmarc (.cmi), 
    affichage des différences, hormis les coordonnées
    """
    raise NotImplementedError

import gc

def find_names(obj):
    """Renvoit les noms de variable de l'objet 'obj' passé en paramètre"""
    frame = sys._getframe()
    for frame in iter(lambda: frame.f_back, None):
        frame.f_locals
    result = []
    for referrer in gc.get_referrers(obj):
        if isinstance(referrer, dict):
            for k, v in referrer.iteritems():
                if v is obj:
                    result.append(k)
    return result

def structuredGrid(nx,ny,nz, delta=-1):
    """
    retourne un tableau np.ndarray de shape(nx,ny,nz,3)
    
    produit tensoriel de X=(x_1,x_2,...x_nx), Y=(y_1,y_2,...y_ny)
    et Z=(z_1,z_2,...z_nz), aléatoirement espacés, mais triés.
    ou bien régulièrement espacés
    Contrainte : nx>1, ny>1, nz>1
    """
#    nx,ny,nz = shape
    if delta == -1 : 
        X = np.random.random(nx)
        Y = np.random.random(ny)
        Z = np.random.random(nz)
        X.sort()
        Y.sort()
        Z.sort()
    else :#TODO: delta=dx,dy,dz
        dx, dy = 1.0/(nx-1), 1.0/(ny-1)
        X = np.arange(0.0, 1.0+dx, dx)
        Y = np.arange(0.0, 1.0+dy, dy)
    if nz == 1 :
        Z = np.asarray([0.0])
    else : 
        dz = 1.0/(nz-1)
        Z = np.arange(0.0, 1.0+dz, dz)
    G = np.ndarray((nx,ny,nz,3), dtype=float)
    for i,x in enumerate(X) :
        for j,y in enumerate(Y) :
            for k,z in enumerate(Z) :
                G[i,j,k,:]=[x,y,z]
    return G

def my2dPlot(XYs, equal=False, cosmetic=('b-','r-', 'b+', 'g*', 'r*','r^')):
    '''
    Tracer des points en 2d
    XYs est une liste de (6 maximum) tableaux de points chaque tableau est de shape (n,2) ou (n,3),
    comporte n points en 2d ou 3d. seules les deux premieres coord de chaque point sont prises en compte.  
    '''
    
    import math
    from matplotlib import pyplot
#    matplotlib.use('MacOSX')
    xmin = ymin =  10^10
    xmax = ymax = -10^10
    colors = cosmetic
    for k, xy in enumerate(XYs) :
#        print xy.shape
        pyplot.plot(xy[:,0], xy[:,1], colors[k])
        xmin = min(xmin, np.min(xy[:,0]))
        ymin = min(ymin, np.min(xy[:,1]))
        xmax = max(xmax, np.max(xy[:,0]))
        ymax = max(ymax, np.max(xy[:,1]))
    w, h = xmax-xmin, ymax-ymin
    dizaine = int(math.log10(w))#10^dizaine <= w < 10^(1+dizaine)
    ax = pyplot.axes()
    if equal : ax.set_aspect('equal')
    ax.set_xlim(xmin-w/10, xmin+w+w/10)
    ax.set_ylim(ymin-h/10, ymin+h+h/10)
    ax.grid(which='major', axis='x', linewidth=0.75, linestyle='-', color='0.75')
    ax.grid(which='minor', axis='x', linewidth=0.25, linestyle='-', color='0.75')
    ax.grid(which='major', axis='y', linewidth=0.75, linestyle='-', color='0.75')
    ax.grid(which='minor', axis='y', linewidth=0.25, linestyle='-', color='0.75')
#    ax.xaxis.set_major_locator(pyplot.MultipleLocator(10**dizaine))
#    ax.xaxis.set_minor_locator(pyplot.MultipleLocator(10**(dizaine-1)))
#    ax.yaxis.set_major_locator(pyplot.MultipleLocator(10**dizaine))
#    ax.yaxis.set_minor_locator(pyplot.MultipleLocator(10**(dizaine-1)))

    pyplot.show()

def load(filinname):
    try : 
        filin = open(filinname,'r')
        return pickle.load(filin) 
    except IOError : 
        print>>sys.stderr, "Impossible d'ouvrir le fichier dump %s, pas de lecture."%filinname
#    raise PyglideIOError("Impossible d'ouvrir le fichier dump %s, pas de lecture."%filinname)
# def doc(object):
#     for key, value in object.__dict__.iteritems() :
#         if type(value).__name__ in ('function',) : 
#             print "%s:\n"%(key)+len(key)*'='+ "\n%s"%value.__doc__
def prettystring(dico, level=0):
    if not isinstance(dico, dict) : 
        return str(dico)
    indic = '\n==>' if level==0 else '\n    *'
    strings = []
    for key, value in dico.iteritems() :
        strings.append(indic+" %s : %s"%(key,prettystring(value, level=level+1)))
    return ''.join(strings) 
#         print "%s:\n"%(key)+len(key)*'='+ "\n%s"%str(value)
def prettyprint(dico):
    '''Pour afficher des dictionnaires (éventuellement imbriqués) de manière lisible'''
    print prettystring(dico)

if __name__=="__main__":
#     help()
    points = np.asarray((np.arange(5), np.arange(5,10))).reshape((-1,2))
    print points
    print hardScale(points, (2,0), centre=None, translate=False)
    
    dico = dict(locals())
    for key, value in  dico.iteritems():
        if type(value).__name__ in ('function',) : 
            print "%s:\n"%(key)+len(key)*'='+ "\n%s"%value.__doc__
    import matplotlib
    app = QtGui.QApplication(sys.argv)
    print phrase()
    from config import VALIDATION_DIR
    filename = Path(VALIDATION_DIR, 'dxf', 'cot_d03.DXF')
    lines = open(filename).readlines()
    print findAll('LINE', lines)
    p1 = QtCore.QPointF(1,0)
    p2 = QtCore.QPointF(1,1)
    print dist2(p1, p2)
    p1 = [1,0.]
    p2 = [1,1]
    print dist2(p1, p2)
    #### segmentPlusProche(points, P):####
    P = np.asarray([(3.0, -1.0),(2,-1),(0,0),(2,2), (3,2), (5,2), (10,0)])
    print P
    P = np.asarray([(0, -1.0),(-1,0),(1,0),(0,1)])
    print 'barycentre', baryCentre(P)
    exit()
    a = np.asarray([2.0, 0.0])
    i = segmentPlusProche(P, a)
    a.shape = (1,2)
    my2dPlot((P,a), equal=True, cosmetic=('b-', 'bo', 'r-', 'ro', 'g*', 'r*', 'r^'))
    print P[i], P[i+1]
#    exit()
#    print rreplace("kze,$$$$$l$$$izefh","$$","$")
#    datadir = '/Users/puiseux/Documents/workspace/pyglide/test'
#    filename = Path(datadir, 'simple.OUT')
#    lines =  open(filename).readlines()
#    print find('yyy', lines)
#    print findAll('xxxxxx', lines)
#    tag = 'TIME STEP'
#    deb = _findRel(tag, lines)
#    print "relatif", _findRel(tag, lines, 100)
#    print "absolu",find(tag, lines, 4000)
#    print findAll(tag, lines)
#    points = np.arange(21, dtype='float')
#    ppp = points
#    points.shape = (7,3)
#    print points
#    print find_names(points)
#    print structuredGrid(2,3,1)
#    print 'doublons : '
#    nbp = 10
#    a = np.random.randint(0, 3, 3*nbp).reshape((-1,3))
##    voisinages = [range(i,nbp) for i in range(nbp)]#[[0,...n-1],[1,...n-1],[2,...n-1],...[n-1]]
#    voisinages = [[i]+range(nbp) for i in range(nbp)]#[[0,0,...n-1],[1,0,...n-1],[2,0,...n-1],...[n-1, 0,...n-1]]
##    print a
#    for k0,k1 in doublons(a, voisinages=voisinages, sym=False) : 
#        print k0, k1
#        print a[k0],a[k1]
#    
    
#    print encombrement(points)
    
#    a = 100.0
#    points = [(-a,-a),(a,-a),(a,a),(-a,a)]
#    print encombrement(points)
#    tags = ['zozo','TIME STEP','COMPONENT COORDINATE SYSTEMS','AERODYNAMIC DATA FOR PATCH','UNIT NORMAL VECTORS AND PANEL AREAS FOR PATCH']
#    
#    print mark(tags, lines)????
    sys.exit(app.exec_())

